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Abstract. — We consider a numerical scheme for Hamilton-Jacobi equations 
based on a direct discretization of the Lax-Oleinik semi-group. We prove 
that this method is convergent with respect to the time and space stepsizes 
provided the solution is Lipschitz, and give an error estimate. Moreover, wc 
prove that the numerical scheme is a geometric integrator satisfying a discrete 
weak-KAM theorem which allows to control its long time behavior. Taking 
advantage of a fast algorithm for computing min-plus convolutions based on 
the decomposition of the function into concave and convex parts, we show that 
the numerical scheme can be implemented in a very efficient way. 



and where uq is a given global Lipschitz function. 

We will mainly consider the case where H is separable, in the sense that we 
can write H(t,x,p) = K(p) + V(t,x), for some convex function K and some 
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1. Introduction 



We consider Hamilton-Jacobi equations of the form 

(1) d t u + H(t, x, Vm) = 0, u(0,x) = u 
where H(t, x, v) is a Hamiltonian function 

(2) H : R x 1" x R™ -> R, 
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smooth and bounded function V . The typical cases of study we have in mind 
are the so called mechanical Hamiltonians, of the form 



where P G 1™ is a given vector, \v\ = uf H h for v = (v%, ■ ■ ■ , v n ) G M n , 

and where V(t, x) is a suitably smooth and bounded function. 

Since the pioneering works of Crandall and Lions [CL84] and Souganidis 
[Sou85], the study of numerical schemes for the Hamilton-Jacobi equation (1) 
has known many recent progresses, see for instance [LS95, LT01, Abg96, 
JX98, BJ05] and the references therein, and more specifically [OS88, OS91, 
JPOO] for the popular (weighted) essentially nonoscillatory (ENO and WENO) 
methods which are now widely commonly used in many application fields. 

Following a different approach, more in the spirit of [FF02, Ror06] (or 
[AGL08] in an optimal control setting), the main aim of this paper is to show 
how a direct discretization of the Lax-Oleinik representation of the viscosity 
solution of (1) allows to define a new fast algorithm for computing u(t, x) pos- 
sessing strong geometrical properties allowing to control its long time behavior 
and obtain error estimates when the solution is Lipschitz. 

Let us recall, see [Lio82, Fat05], that under some assumptions on H 
(smoothness, uniform superlinearity and strict convexity over the fibers, see 
Section 2 below), we can write 



where the infimum is taken over over all absolutely continuous curves 7 : 
[0, t] —?• W 1 such that j(t) = x, and where L(t, x, v) is the Lagrangian asso- 
ciated with H. The idea of this paper is to discretize directly (4) on a space 
time grid, by replacing the set of curves 7 by the set of piecewise linear curves 
across the space grid points. 

We first prove that such an approximation is convergent with respect to the 
size of the space and time stepsizes, and under an anti-CFL condition (namely 
that the ratio between the space and time stepsize should be small). We give 
an error estimate under the assumption that uo is Lipschitz. 

Moreover, this numerical integrator turns out to be a geometric integrator 
(see for instance [HLW06, LR04]) in the sense that is respects the long time 
behavior of the exact solution u(t, x). Let us recall that in the case of periodic 
Hamiltonians (both in time and space variables), the weak-KAM theorem (see 
[Fat05, CISMOO]) shows the existence of a constant H such that 



(3) 




(4) 



u(t, x) 




1 



u(t, x) — > H when t — > +00. 
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Here, using a discrete weak-KAM theorem, see [BB07, Zavl2], we prove that 
the numerical scheme possesses the same long time property, with a constant 
that is close to the exact constant H. 

Finally, we show that in the separable case mainly considered in this pa- 
per (see (5) below), the discrete version of (4) is a min-plus convolution that 
can be approximated using a fast algorithm with 0{N) operations in many 
situations if N is the number of grid points. This algorithm uses the decom- 
position of u into concave and convex parts. Moreover, it easily extends to 
any space dimension n using a splitting strategy, when the kinetic part of the 
Hamiltonian is separable - see Remark 2.8 - which includes the case (3). 

We then conclude by numerical simulations in dimension 1 to illustrate 
the good behavior of our algorithm, as well as its very low cost in general 
situations. 

The paper is divided into three parts: in a first part (Section 2) we give a 
convergence result over a finite time interval of the form [0, T] where T is fixed. 
In a second part (Section 3), we consider the case where the Hamiltonian is 
periodic in time t and x. In this case, we can derive explicitly the dependence 
in T in the error estimates, and prove a weak-KAM theorem for the numerical 
scheme which gives informations concerning the long time behavior of the 
scheme. In the third part (Section 4), we describe the implementation of the 
method based on a fast algorithm to compute min-plus convolutions. We 
conclude this part by showing numerical simulations. 

Acknowledgement. — This work owes a lot to Vincent Calvez, who put 
the authors in touch and took part to preliminary discussions. It is a great 
pleasure to thank him a lot. We also would like to thank Vinh Nguyen for 
careful reading through previous versions of the paper. Finally, the last author 
would like to thank Antonio Siconolfi for bringing him to this subject. 

2. Description of the scheme and convergence results 
2.1. Hypotheses. — We consider a Hamiltonian H(t,x,p) of the form 



defined onix 1" x W 1 . With this Hamiltonian we can associate by Legendre 
transform the Lagrangian 



(5) 





and we calculate that in our case 



L(t,x,v) = K*(v) - V(t,x) 
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where K*{v) is the Legendre transform of K. For instance in the special case 
(3) we have 

1 2 

L(t,x,v) = -\v\ — P-v — V(t,x). 

We make the following assumptions on K and V: 

(i) The function K* E C 2 (M™) is uniformly strictly convex in the sense that 
there exists a constant c > such that for all Y E W 1 , and for all v E M n , 

(6) ^(v)(Y,Y)>c\Y\ 2 . 

(ii) The function V(t, x) E C 2 (M x W 1 ) is such that there exists a constant B 
such that for j + q ^ 2, and all (t, i) eRx R n , 

(7) 

where | • | denote the norm of differential operators acting onlx M. n . 

Note that the bound (7) is straightforward under the supplementary as- 
sumption that V(t,x,v) is periodic in (t,x), the case studied in the next sec- 
tion. 

Remark 2.1. — The previous hypotheses imply that the Hamiltonian H and 
the Lagrangian L are C 2 , convex and superlinear in respectively p and v: 

(8) Vfc>0, Vt>0, 3,40) <oo, L(t,x,v) ^ k\v\ - A(k). 

Under these assumptions, the viscosity solution of (1) can be represented 
by the formula: for all t, 6 > 0, 
(9) 

rt+S 

Vi £ 1", u(t + 5,x) :=T t s u(x) = inf u(t,j(t))+ L(s,j(s),j(s))ds, 

-y(t+S)=x J t 

where the infimum is taken on all absolutely continuous curves 7 : (t, t + S) — > 
W 1 verifying j(t + 5) = x, see [Lio82, Fat05]. Moreover, the infimum is 
achieved on a curve jf x (s) that is C 2 and satisfies the Euler-Lagrange equation 

(10) ~~( S , 7 (s),7(s)) = ^(*,7(«).7«)- 

as ov ox 

The notation T/ defines the Lax-Oleinik semi-group. In particular, we have 
Tg^g o T/ = Tf +(T for non negative 5 and a. With these assumptions, we have 
the following Proposition. 

Proposition 2.2. — For all T > 0, and for all R > 0, there exists M(R,T) 
such that for all x, y £ M. n satisfying \x — y\ ^ R and for all t E K, then every 
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solution of the Euler- Lagrange equation (10) minimizing the action 

rt+T 

(11) J L{s, 7 (s),j(s))ds, 

with fixed endpoints *y(t) = x and j(t + T) = y, satisfies \^{s)\ ^ M(R, T) for 
all s G [t,t + T]. 

Proof. — The Euler-Lagrange equation is written 

Using the uniform strict convexity of K* and the fact that d x V is uniformly 
bounded, there exists a constant C depending only on T and K, such that 

(12) Vs£[t,t + T] | 7 (s)| 2 ^C. 
This implies that for all s E [t,t + T], 

(13) | 7 (s) - j(t)\ < [ t+T \j(s)\ ds TVC. 



Now as 7 minimizes the action between t and t + T, comparing with the trivial 
curve 1 1->- x + t(y — x)/T from x to y, we get 

L(sMs),i(s))ds ^ J K*( y —)-V(s,x + -{y-x))ds 

^ TD(R/T) + TB, 

where D(M) = sup^^^- |ff*(u)|, and B is given by (7). By superlinearity, we 
deduce that 

rt+T 

^ \i(s)\ds ^TD(R/T) + TA(1) + TB, 
where A(l) is given by (8). Using (13), we thus obtain 

T| 7 (t)| ^ J + \j(s)\ds + T 2 VC ^ TD(R/T) + TA(l) +TB + T 2 Vc. 

This shows that |7(t)| is bounded, and hence using again (13) that (7(5) | is 
bounded for all s, with a constant depending only on T, R and the constants 
appearing in (6) and (7). □ 

Remark 2.3. — The previous lemma is one of the main keys in the proof of 
the convergence of our schemes. Here, it is established thanks to the particular 
form of H, but it can be noted that it remains valid under other technical 
assumptions (for example if H is autonomous and Tonelli as established in 
[Fat05, FM07], or if it is Tonelli and periodic both in the space and in the 
time variable as proven in [Mat91, CISMOO, Itu96]). Actually, in these 
cases, it can be established that M(R,T) only depends on the ratio R/T. 
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We will come back on these matters in Section 3 and give a proof of this 
result in the Appendix. Therefore, this section and the next would still be 
valid for general Hamiltonians chosen in these classes, however the convolution 
techniques of section 4 would fail. 

Finally, a clear consequence of Equation (12) is that the Euler-Lagrange 
flow of L is complete. 

2.2. An approximate semi-group. — For a given e > we define the 
e-grid G e = eTU 1 endowed with the metric induced by the euclidian metric on 
W 1 . For a given continuous function u, we define u\g e '■ G £ — > R its restriction 
to the grid G e . Given t, r > 0, let us define c^ £ : G £ h- > E as follows : 

(14) V(x,y)GGi cl £ (x,y) = j^\^ x + ( s -t)V^ V -^j&s. 

Let us introduce the following discrete Lax-Oleinik semi-group: if u : G e — > 
K is any function, we set 

Vx e G e , T t ' u(x) = inf u(y) + c£ Jy,x). 

A good setting to apply this semi-group is the one of functions u with 
linear growth, which means that the quantity (u{y) — u(x))/(l + ||y — x||) 
is uniformly bounded. Note that for such a function u, for a given x, the 
hypotheses made on L ensure the existence of a minimizing y G G £ attaining 
the previous infimum. Indeed, the cost cj £ inherits from L a superlinearity 
property which implies that the function y \-t u(y) + cl £ (y, x) goes to oo at oo. 
Moreover, the set of functions with at most linear growth is invariant by our 
semi-group. For more details, we refer the reader to the appendix of [Zavl2]. 
Note that Lipschitz functions, which we will only consider in the following, 
have linear growth. 

For a given integer N, we define 

(15) T t fu = T[ N _ 1>£ o.-.o Tl £ o T[ (h£ u 

the composition of N times the discrete semi- group T[ £ , where for all i = 
1,...,N- 1, U = t + ir. 

Remark 2.4- — The following monotonicity property holds true: if u(x) ^ 
v(x) for all x E M n , we easily observe that T^Ju ^ T^Jv. 

Proposition 2.5. — Let u : W 1 \— > M and N ^ 1 an integer, t E IR and r > 0. 
Then 

(16) (T t NT u)\ G£ ^T t f(u\ Ge ). 
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Proof. — Let x G G £ . With this point, we can associate N + 1 points Xi, 
i = 0, . . . , N such that xm = %, and such that for all i = 0, . . . , N — 1, 

T ^, £ ° T t>(^+i) = T£u(xi) + c[ ii£ (xi,Xi + i). 
With this notation, we can verify that 

N-l 

T t N e T ( u \G £ )(x) = u(x ) + CtiA Xi > Xi +^ 



N-l rt-.n 



2=0 

Xi-\-~l Xi Xi-\-\ Xi 



= u(x ) + ^2 / L is,Xi + (s -ti)^ 1 ^ ^^+J l -\ds 

pNt 

= u(x ) + / L(s,7 ti£)T (s),7 t)E)r (s))ds, 

J 

where 7t ie , r (s) is the continuous piecewise linear curve defined by 

(17) j t£T (s) = x i + (s-t i ) a ^ — — , for se[ti,t i+1 ). 

T 

The result then easily follows by definition of T^ T u. □ 

As we will see now, the reverse inequality (16) is true up to a small error 
term coming from the time and space discretization. This is stated in the 
following convergence result: 

Theorem 2.6. — Let T > 0, e$, tq and ho > and u : W n — > R a bounded 
Lipschitz function. There exists a constant M such that for all e > and 
t > such that e < Eq, t < tq and 

(18) - < ho, 

T 

for all N satisfying Nt ^ T , 

(19) (T t N -u)\ Ge -T t N /(u\ Ge ) ^M( £ - + r). 

oo T 

Proof. — Let us denote by jt ■ [t,t + T] — > M n a minimizer of (9). Recall that 
the curve 7t(s) is C 2 . Let us set y := 7*(i) and x := jt(t + T). We have 

rt+T 

T?u(x) = u(y) + J L(s,7t(s),7 t (s))ds. 
By superlinearity (8), this implies that 

rt+T 

\jt(s)\ds ^ TA(1) + \T?u(x) - u{y)\. 
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Comparing with the trivial curve 7 = x in the definition of the Lax-Oleinik 
semi- group, we have that 

(20) T^u(x)^\u\ oo + T(B + K*(0)) 

where B is the constant in (7). Moreover, since L is bounded below 
(L(t,x,v) ^ b for some constant b, for all (t,x,v)), clearly, the action of any 
curve defined for a time T is greater than Tb which implies immediately that 

(21) Tlu(x)^-\u\ oo + Tb. 

Hence, there exists a constant B\ depending only on L and \u\ such that 

rt+T 

(22) \ x - y \ = \ lt (t + T)- lt {t)\^ J |7 t («)|da<Bi(l + r). 

Remarking that 74 is also a minimizer of the action (11) under the constraint 
j(t) = y and 7(2 + T) = x, we can apply Proposition 2.2 which shows that 
there exists a constant M\ = M{B\(\ + T),T) depending only on T, L and 
I ii I such that 



■ 00 



(23) VsG[M + T], \j t (s)\^Mx. 

Assume now that N is an integer such that Nt ^ T. For all % = 0, . . . , N 
we define 

1 , 

-7t(*» 

where for i = 0, . . . , N, ti = t + it and where the function |_ - J is the floor 
function, coordinate by coordinate. With these points, we associate the con- 
tinuous piecewise linear path 7$ eT defined as in formula (17). Notice however 
that the points X{ are no longer the same. By definition of the points Xj, we 
have 



rn+i 

^2ey/n + \jt(s)\ds sC 2ey/n + rM 1 . 

Ju 



(24) Wie[0,N], \x i -j t (t i )\ = \j ttEjT (t i )-jt{ti)\<£Vn. 
Now, using the bound (23), we have for all i = 0, . . . , N — 1, 

But this inequality implies that for alH = 0, . . . , N — 1, 

Vs G [ti,t i+ i], |7t, e , T (s) - Xi\ < 26^ + rMi, 
while |7t(s) — 7t(£i)| ^ tMj upon using (23). Hence we get 

(25) Vs € Mi+iL |7t,e,r(s)-7t(*)l < 3eVn + 2TMi. 
Moreover, we have for s,a & [t< , 5 
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upon using (12). Hence for s G [£i,ti+i], we have 

\jt(U+i) - lt{U) - rj t (s)\ < / |7*(cr) - j t {s)\da < t 2 C 



and hence for all s S [t^, 

7t(*i+l) -7*(<i) 



7H* 



^ rC. 



Using (24), we obtain easily that for alH = 0, . . . , N — 1, 



(26) 



2s 



n. 



Note that using (18) and (23), the previous equation implies that for all i 
Q,...,N-1, 



(27) 



VsG [U,ti + i], \% E>T (s)\KM 2 



for some constant M2 = tqC + h^^fn + M\ independent of e and r. 
Now by definition of Jt,e,r, we have 



r-t+Nr 

(28) / L(s,7 t (s),7t(s))ds - ^ c[ ii£ (xi,x i+ i) 
^* i=0 



Vr 



£(s,7t( s )>7t( s )) ds 



A'r 



^(S,7t,£,r(s),7t,£,r(s))ds 



8=0 



Hs,lt(s),jt(s)) - L(s,Jt,e,T(s), A ft,e,r(s)) 



ds. 



Using (7), the fact that K* is C 2 and the bounds (23) and (27), there exists a 
constant M3, depending on L, M\ and M2, such that the previous error term 
is bounded by 

AT-i , fi+1 

M 3 £ / " 7t, £ ,r(s)| + |7t00 " 7t >e>T (a)|)ds- 

i=0 

Using (25) and (26), this shows that there exists a constant M4 independent 
on e and r, such that 



(29) 



t+Nr N-l 

i=0 



«S M 4 (- + r) 



where we used the fact that e ^ tqe/t. 
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Now by definition of T( £ , we have using (24) and the Lipschitz nature of 

u, 

N-l 

T t N 6 T ( u \G,){x) «S u(x )+J2 C U,e( X ^ X i+^ 

^ (T t NT u)\ Ge (x) + |«(x ) - «(7t(0)) | + M^- + r) 

(30) < { T t ^u)\ Gs { X ) + M{ £ - + r) 

for some constant M independent on e and r. This inequality and (16) show 
the final estimate (19). □ 

2.3. Fully discrete semi-group. — In the previous Section, we have seen 
that the Lax-Oleinik semi- group can be approximated using the cost (14) 
defined on the grid. To compute this cost, a quadrature rule in time has to 
be used. In this subsection, we prove how the Euler approximation of this 
integral yields a convergent scheme which still satisfies a weak-KAM theorem 
similar to Proposition 3.1 under suitable periodicity assumptions. 
For a given e > and r > 0, we define the following cost function: 

(31) \/(x,y)eG 2 £ , KL(y,x)=TL(t,x, X ~ !J 



T 

and the associated fully discrete Lax-Oleinik semi-group 

Vx G G £ , T t T £ u{x)= inf u{y) + K T t (y,x), 

if u : G £ — > M is a function. Using the explicit expression of L, we can rewrite 
this fully-discrete semi-group as 



(32) Vx G G e , T t ] £ u{x) = inf \u{y) + tK* \ - rV(t,x). 

involving the (min,plus)-convolution of u and K* . 

Remark 2. 7. — We can interpret this scheme as a discretization of the split- 
ting scheme (see for instance [JKR01]) with time step r based on the decom- 
position 

d t u(t,x) + K(Vu(t,x)) = Q, and d t u(t, x) + V(t, x) = 0, 

where the first part is integrated using the method described in the previous 
section. 

Remark 2.8. — In dimension n 1, if we assume that for p = {pi, . . . ,p n ) G 
M n , K(p) = K\(pi) + • • • + K n {p n ) with convex Hamiltonian functions K* , 
i = 1, . . . , n, satisfying all the hypotheses (i), (ii) on M, then we immediately 
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see that for a given function u(x) = u(x\, . . . , x n ), with x = (xi, . ■ ■ , x n ) £ M n , 
we have 



inf u (y) + rK*(^-) 



inf 



inf tK, 



T 



+ 



inf tK*(— — — ) + u{y 1 ,...,y n ) 



where we have decomposed G e = G\ x • • • x G™ and where 



Vi G [l,n], T??u{x) = inf . T K*{^^) + 



T 



u{x\ , . . . , Xi—\ , y^, Xj+i , . . . x r 



This formula is essentially due to the fact that the Hamiltonians Ki commute, 
i.e. satisfy {Ki,Kj} = for (i, j) E {1, . . . , re} 2 which ensures that the flows 
of dtu = Ki(Vu), i = 1, . . . ,n commute. In this case, this allows to reduce the 
computation of the minimum over the n dimensional grid G £ to n minimization 
problems over the 1 dimensional grids G\. 



For a given integer N, we define - compare (15) 

T t N £ T u = V N _ 1 , £ o---oT t l £ oT t r , £ u 



where ti = t + it. Note that with these notations, an estimate of the form 
(16) is no longer valid. However, we have the following convergence result: 



Theorem 2.9. — Let T > 0, Eq, tq and h$ > and u : M. n — > R a bounded 
Lipschitz function. There exists a constant M such that for all t > 0, and all 
e > and r > such that e < eq, t < tq and the bound e/t < ho are satisfied, 
then for all N verifying Nr ^ T , 



(T t NT u)\ Ge -T^(u\ Ge ) 



Nt, 



sC M(-+t). 

oo T 
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Proof. — The proof is exactly the same as the proof of Theorem 2.6 until 
Equation (28) that has to be replaced by the following estimate: 

t+Nr N-l 

L(s,j t (s),jt(s))ds - Kl £ (xi,Xi +1 

i=0 

JV-1 



L[s,-y t (s),jt(s)) - L(s,7 Mi7 -(s),7 tj£jT (s)) 



ds 



i=0 



i(s,7*,£,r(s),7t,£,r(s)) - i(ii,7t,e 1 r(ii),7t,e,r(ii)) 



ds. 



In the right-hand side of this inequality, the first term is bounded by M4 (~ +r) 
see (29). To bound the second term, we observe first that for s 6 [ij, ii+i], the 
derivative 7t )£)T (s) = (xi+i —Xij/r does not depend on s. Hence using (7) and 
(27) the function 

[ti,ii + i] 9s4 L(s,7 tj£)r (s),7 t)£)T (s)) 

is C 1 with uniformly bounded derivative. Thus we obtain that there exists a 
constant M5 such that 

£(s,7i,e,'7-( s )>7t,e,'r( s )) - L(t i} 7i, e , r (ti), 7t, £ ,r(^)) < M 5 (s - ij). 

This proves that (compare (29)) 



(33) 



i=0 



<Af 6 (- + r), 



for some constant Mq independent of e and r. As in (30), we obtain that 

T t N /(u\ Ge )(x) < (If T «)|o.(x) + M(J + r), 

for some constant M independent on e and r. 

To prove the reverse inequality, let us fix x 6 G e . We consider a sequence 
i = 0, . . . , iV with un = x and 

JV-1 

(34) r t N /(u\ Ge )(x) = u(y ) + ^ K l >£ ( yi ,y i+1 ), 

8=0 

and we define the curve 

Vt,eA s ) =yi + (s-t i )^ : — — , for se[ti,t i+1 ). 

T 

Note that in a similar manner to what we did to prove the inequalities 20 and 
21, using the fact that u is bounded, and comparing with the trivial sequence 
made of a constant point (with (7)) on the one hand, and the fact that L, 
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hence the kJ , are bounded below on the second hand show that there exists 
a constant D\ such that 

|72 T («k)L<0i- 

By superlinearity of L (and of the kJ. £ ) and using again the fact that u is 
bounded, we thus see, as in 22 that there exists a constant D 2 such that for 
alH = 0, .. . ,N- 1, 

' Vi+i ~ Vi 



^D 2 , 

T 

which in turn implies that 

Vs £ [ti,t i+ i], \r]t, £ ,T(s) ~ r] t ,e,r(ti)\ ^ tD 2 . 

As the derivative of r] tt£tT (s) with respect to s is uniformly bounded by D 2 and 
constant on the time intervals [t^tj+i], and as L is C 1 with uniformly bounded 
derivative on M. x M. n x B(0, D 2 ), we obtain 



(35) 



iV-l 



i=0 



Nt 



^2 ^t^eiVi+l^i) ~ / L(s,rft, £ , T (s),r) t ,eA s )) ds 



for some constant D3. Using the definition of the exact semi-group, we thus 
have 



Nt 



( T t NTu )\G £ ( x ) < ^(^,e,r(*jv)) + / L(a,rft, B , T (a),Tft^ T (s))ds 

Jo 

< ^ t («|gJ(x)+tL» 3 

upon using (34) and (35). This proves the result. □ 



3. Long time behavior in the periodic case 

We will now make the supplementary assumption that the potential function 
V(t, x) is periodic, namely 

(Hi) The function V is Z X Z n -periodic, in the sense that 

V(f,i)eKxl", V(m, M) eZx Z n , V(t, x) = V(t + m, x + M). 

Moreover, V £ C 2 (l x R n ). 
Note that under this assumption, the estimate (7) is automatically satisfied. 
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3.1. Weak-KAM theorem and a priori compactness. — In this peri- 
odic case, the weak-KAM theorem allows to study the long time behavior of 
the solution of (1) defined by the Lax-Oleinik semi-group: 

Proposition 3.1. — Assume that the hypothesis (i) and (Hi) are satisfied. 
Then there exists a unique constant H such that there exists aZx Z n -periodic 
continuous function u* : M x R™ — > E verifying for all t > 0, 

T t 1 u*(t,-) = H + u*(t,-). 

Moreover, for any uniformly bounded u : M. n —> M., there exists a constant C u 
such that we have 

Vt>0, |T *u - tH]^ < C u . 

Proof. — The existence of H and of u* is exactly the content of the weak- 
KAM theorem (see [Fat97, Fat05] for the autonomous case, and [CISMOO] 
for the time periodic case). The second assertion is a consequence of the fact 
that 

\l$u-7$v\__ ^ \u-v\ 



1 oo 1 1 oo 



for all continuous bounded functions u and v on l n , where | • |oo denotes the 
L°° norm on W 1 . □ 

The goal of this section is to prove that a similar result holds for the nu- 
merical scheme described in the previous section, under the assumptions (i) 
and (Hi). 

In order to study the long time behavior of the method in this case, we 
first give an a priori compactness result which refines the estimates given 
in Proposition 2.2. The following proposition is mainly due to Mather (see 
[Mat91] for the case of time periodic Lagrangians, or [Itu96, Lemma 7 and 
Corollary 8] for space periodic Lagrangians). 

Proposition 3.2. — Assume that H satisfies the hypotheses (i) and (Hi). 

For all r > 0, there exists a constant V such that for any minimizer of the 
Lagrangian action 7 : [a, b] — > W 1 with b — a ^ 1 and (7(6) — 7(a) \/(b — a) ^ T 
then we have 

ViG [a, 6], |7(i)|^T'. 

In other words, the constant M(R, T) of Proposition 2.2 can be chosen to 
be an increasing function of R/T. 

For the sake of completeness, we will give in appendix a proof of this propo- 
sition. Note that most of the proof - essentially taken from [Mat91] - does not 
require the Hamiltonian to be periodic in space. In [Itu96] a similar result 
is proven which requires the Lagrangian to be periodic in space, but not any 
more in time. 
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3.2. Convergence estimates in the periodic case. — Using the previ- 
ous proposition, we can compute explicitly the time dependence in the error 
estimates of Theorems 2.6 and 2.9 in the periodic case. 

Theorem 3.3. — Let To > I, eq, tq and ho > and u : M n — > R a bounded 
Lipschitz function. There exists a constant M such that for all e > and 
r > such that e < Eq, t < tq and 

(36) -<ho, 

T 

for all N satisfying Nt To, 
and in similar way, 



{T^u)\ Ge - T t N /(u\ Ge ) < MNr( £ - + r), 

oo T 



(37) (T t NT u)\ Ge -% N /(u\ Ge ) < MNr( £ - + r). 

oo T 

Proof. — In the proof of Theorem 2.6, equation (22) then gives using Propo- 
sition 3.2 (with T To > 1) that the constant Mi defined in (23) does not 
depend on T = Nt and depend in fact only on To. It then follows that M2 
and M3 also are independent of T = Nt, while M4 is proportional to the time 
of integration, that is Nr. In a similar way, in Theorem 2.9, the additional 
error made is proportional to the time of integration, hence giving the second 
part of Theorem 3.3. □ 

3.3. Discrete weak-KAM theorem and effective Hamiltonian. — 

Recall that the function u(t, x) is defined on T 1 xT n = (R/Z) x (R n /Z n ) . For 
convenience, we will only treat the cases of rational time and space discretiza- 
tions: We set 

1 1- 

and in the sequel, we will only consider stepsizes (e, r) E A. We then will 
denote by p the canonical projection from W 1 to T n , and by G e = G £ /Z n the 
quotiented grid, where G e is the grid above, defined on M n . 

Finally, we define two new cost functions: For (e, r) G A and t > 0, 

V (x, y) G (G £ ) 2 , c$ ' (x, y) = inf _ cL(x, y), 



A={(-,-)|(M)€N*xN*}, 



p(x)=x 

p(y)=y 

where cj £ (x,y) is defined in (14), and similarly 

V(x,y) G (G e ) 2 , Kj (x,y)= inf Kj e (x,y), 

p(x)=x 

p{v)=v 

where Kj e (x,y) is the fully discrete cost function defined in (31). 
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We then define the following semi-groups: Let u : T n — > E and (e, r) G A, 
then 

(38) Vi G G £ , fl £ u(x) = inf u(y) + ^ £ {y , x) , 

and the fully discrete semi-group 

Vx G G e , T t T £ u(x) = inf tt(y) + K[ £ (y,x). 

As in the previous section, we define 

T ^ T = T t N -ue ° • • • ° T f T llE ° T r o ,e 

and 

% N e T = Vn-^s ° • • • 7^, £ 7^ >e > 

where ij = t + ir. 

We leave to the reader the verification that if u : T n — > IR is a function 
and if u : M. n — > M is its lift (which is then Z n periodic), both functions T^Ju 
and Tt N £ T u are Z n periodic and that the functions they respectively canonically 

induce on T n are Tj^Ju and T t N £ T u. It comes from the fact that two infimums 
commute. Hence the previous convergence result Theorem 3.3 can be read 
equivalently on T n or on the space of Z n periodic functions on M. n . 

We can use the discrete weak KAM theorem to better understand the ap- 
proximate semi-groups applied for a period 1 of time and obtain the following 
proposition. 

Proposition 3.4- — For any (e,r) G A, there exist unique constants H £)T 
and % £ ,t such that there exist functions u* T , v * T : G e — > M verifying: 

^0,e u *,T = u *,t + H £tT , 

and 

Moreover, in u is any bounded initial datum on G £ at t = 0, then we have in 
L 



OG 



and 



as N — > +oo. 



Proof. — The first part is just a reformulation of the discrete weak KAM 
theorem (see for example the appendix of [Zavl2] or [BB07]) while the second 
part is - as in the proof of Proposition 3.1 - a direct consequence of the fact 
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that our approximation operators are weakly contracting for the infinity norm 
on bounded functions. □ 

Remark 3. 5. — Note that the previous proposition can be interpreted in the 
(min,plus) framework: Equation (38) is the (min,plus) product of a vector u 
by the matrix c[ e . Then, for r = l/£, I £ N*, Tq £ u = Tq t £ u is obtained by 
successive matrix multiplications with u. Hence there exists a matrix C £ , T 
such that Tq £ u(x) = inf ~ g( ^ u(y) + C £)T (y,x), with C £)T (y,x) < +oo for all 
y, x. The matrix C £)T has a unique eigenvalue H £tT , and u £)T it an eigenvector 
(see [BCOQ92] for' details). 

Let us recall that H is the effective Hamiltonian of if. It is obtained in 
homogenization theory by solving the cell problem ([LPV87]) and is also the 
constant found in the weak KAM theorem (3.1). Using the refined convergence 
result obtained in Theorem 3.3, we can estimate the error between the effective 
Hamiltonian and the discrete effective Hamiltonians defined in Proposition 3.4. 

Theorem 3.6. - - With the notations of Theorem 3.3, let (e,r) £ A be such 
that e ^ eo, r ^ tq and e/r ^ ho, then following inequalities hold: 



\H £iT -H\ ^M(-+r), 



r 



and 



\n e , T -H\ ^ M(- + t), 



r 

where M is any constant coming from 3.3, and where H £)T and% £)T are defined 
in Proposition 3.4-. 

Proof. — We will only prove the first inequality, the second being obtained 
in the same way. Start with a bounded and uniformly Lipschitz continuous 
function u : R n — >■ M. By 3.3, the following inequality holds if Nt ^ To for 
some chosen Tq > 1 : 



(39) (T NT u)\ Ge -T ^(u\ Ge ) 



€ MNt(- +t) 



T 



Dividing by Nt and letting TV go to oo yields that 

\He,r-H\ <^M(-+t). 



□ 



Remark 3.7. — By Proposition 2.5, we always have the following additional 
information: H F T < H. 
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Remark 3.8. — One may wonder what is the behavior of the quantity 
(Tq t u) I g e — T^(u\g e ) as T = Nt — >• oo. The previous results show that it 
has a linear growth, of rate H — H £jT . Comparing with weak KAM solutions 
yields that the second error term is always bounded. However, in some cases 
more can be said. Indeed, in the autonomous case (L independant of t) 
Fathi proved the convergence of the Lax-Oleinik semigroup ([Fat98]), that 
is, for any initial condition u there exists a weak KAM solution u* such that 
(T^ T u) — NtH — > u* uniformly. Moreover, it can be proved that the iterated 
powers of a (min, +) matrix are periodic after a finite time. Therefore, for 
N big enough, the sequence TqJ (u\g b ), is periodic after a certain time. In 
conclusion, in the autonomous case, one can write 

(T^u)\ Ge - T^(u\ Ge ) = Nt(H - H e>T ) + w N , 

where wn is asymptotic to a periodic sequence. 



4. Fast (min,plus)-convolution 

As we have seen in (32), the numerical scheme considered in this paper 
involves the computation of the (min, plus) convolution 



M,{" [y)+TK ^r))- X ^°'- 

In dimension 1, if the grid G £ is discretized by retaining N points only, the 
numerical cost is a priori of order N 2 . As we will see now, we can use a 
fast (min,plus)-convolution algorithm that turns out to have a linear cost (i.e. 
proportional to N) in many situations. 

In order to ease the presentation, we will not deal with functions defined on 
a grid, but on functions defined on finite and closed intervals. 

Let a,b 6l with a < b. We write / : [a, b] — >• K if / is such that 

f(x) < oo if x E [a, b] 
f(x) = oo otherwise. 

For / : [a, b] — > R, we say that / is respectively convex, concave, affine if 
f\[a t b] is respectively convex, concave, affine. Let / : [a, b] — > M and g : [c,d] — > 
R. The (min,plus)-convolution (or convolution in the remaining of the paper) 
of / and g is defined by: Vi G 1, 

(40) f*g(x) = Mf(y) + g(x-y). 



Recall that f * g = g * f ■ As / and g are finite only on an interval, it is easy 
to see that Vx £ [a + c,b + d], 

f*g(x)= inf f(y) + g(x - y), 
ye[a,b] 
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and Vx ^ [a + c, b + d], f * g(x) = oo. 

We will only consider piecewise affine functions and decompose them ac- 
cording to their affine components: there exists ao = a < a,, < • • • < a n = b 
such that 

/= min fi, 

ie{0,...,n-l} 

where fi : [aj, a^+i] — > M. is an affine function. For i = 0, . . . , n — 1, we denote 
by fi = (/(ai+l) - /(oj))/(oi + i - di) the slope of fi or the slope of / on 

[di, Qi+i]. 



4.1. Convolution. — The fast algorithm to compute the convolution (40) 
is based on a decomposition of g (= u) in piecewise convex and concave func- 
tions. As the function / (= K*(-/t)) considered will always be convex (see 
(6)), we thus see that we are led to compute separately the convolution of 
convex by convex functions, and concave by convex functions defined on finite 
intervals. As we will see, each block can be computed at a linear cost. In the 
end, the global cost of the algorithm thus depends on the number of convex 
and concave components on /, a number which might increase in the time 
evolution of the numerical solution of the Hamilton-Jacobi equation. We will 
come back later to this matter, but we emphasize that this procedure can be 
very easily implemented in parallel, each convolution block being calculated 
independently. 

We start with the following result, the proof of which can be found for 
example [BT08]. 

Lemma 1^.1 (convolution of a convex function by an affine function) 

Let f : [a, b] — )■ M be a convex piecewise affine function and g : [c, d] — > M 
be an affine function of slope g' . Then f * g : [a + c, b + d] — > R is a convex 
piecewise affine function defined by 

{f{x — c) + g{c) if a + c^x^a + c, 
f{a) + g(x — a) if a + c<x^a + d, 
f(x-d)+g(d) if a + d<x^b + d, 

where a = min{oj in the decomposition of f ' \ f-^ g'}. 

Figure 1 illustrates this lemma. In the rest of the section, we will use a 
decomposition of such a convolution into three parts '■ f * g = min^ 1 , g c , g 2 ), 
where 

(i) 9 1 = f * 9\[c+a,c+aY, 

(fi) g c = f * g\[c+a,d+oi\; 

(iii) g 2 = f * g\[ d+a ,d+b}- 
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a a b c d a + c a + c a + d b + d 

Figure 1. Convolution of a convex function by an affine function 
and decomposition into three functions. 

In other words, g 1 is composed of the segments of / whose slope are strictly 
less than that of g, g c corresponds to the segment g and g 2 is composed of the 
segments of / whose slope are greater than or equal to that of g. Note that 
g c is also concave. 

A direct consequence of this lemma is the following theorem, stated 
in [LBT01]. A complete proof is presented in [BJT08]. 

Theorem 4-2 (convolution of a convex function by a convex function) 

If f and g are convex and piecewise affine, then f * g is obtained by putting 
end-to-end the different linear pieces of f and g sorted by increasing slopes. 

For sake of completeness, we give below Algorithm 1 for computing the 
(min,plus)-convolution of two convex piecewise affine functions. 



Algorithm 1: Convolution of two convex functions 
Data: / : [0, n] — > M. a convex function with slopes (r^), g : [0, m] — > M. a 

convex function with slopes (pi). 
Result: h = f * g 
begin 

i^0;j^0; h(0) <-/(0) + fl (0); 
while i + j < n + m do 

if i ^ n and (r% < Pj or j = m) then 

j_ h(i + j + 1) 4- h(i + j) +n; % 4- i + 1; 

else 

L h(i + j + 1) «- h(i + j) + Pj ; + 1; 
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We now turn to the case where / is convex and g is concave. We begin with 
the following lemma: 

Lemma 4-3. — Let f : [a, b] — > K be a convex piecewise affine function and 
g : [c, d] — > M be a concave piecewise affine function which decomposition is 
g = min^L 1 gj . Then 

m 

f*g = mmf*gj. 

Proof. — This is a direct consequence of the distributivity of * over the min- 
imum. □ 

Now, consider two functions / : [a, b] — > E, convex, and g : [c, d] — > R, 
concave, with respective decompositions in f, L : [oj,aj + i] — > M, i G {0, n — 1} 
and gj : [cj,Cj+x] — > M, j £ {0,m — 1}. The following lemma, that considers 
two consecutive affine functions of g, leads to an efficient algorithm to compute 
the convolution of a convex function by a concave function. 

Lemma 4-4- — Consider the convolutions f * gj-x and f * gj. Let 
ctj = min{aj in the decomposition of f \ f[ > g'A 

and 

ctj—i = min{aj in the decomposition of f \ f[ > g'j-\}- 

Then 

- Mx ^ Cj +ctj, f * gj(x) ^ / * gj-i(x); 

- Vx ^ cj +aj-i, f *g j - 1 (x) ^ f*9j(x). 

Proof. — First, as / is convex and <7j_i > g'j, we have that ctj—i > ctj. From 
Lemma 4.1, for all x ^ Cj + ay, 

/ * 9j(x) = f{x -Cj)+ g(cj). 

Either x ^ cy_i + ay_i, then / * gj_i(x) = f(x — cy_i) + g(cj-i) and 

/ * 9j(x) - f * 9j-i(x) = f(x - Cj) - f(x - Cj-i) + g(cj) - g(cy_i); 

as x — Cj-i ^ ay-i, then f(x — cy-i) — f(x — Cj) ^ 9j-i ' ( c j ~ c j-i) an d 
f*gj(x) -f*g j _ 1 (x) ^ 0; 

or x > Cj_i + ctj-x, then / * gj^\(x) = /(ay_i) + g(x — ay_i); as cy_i < 
x — ay_i ^ x — ctj ^ cy, p(cj) — — dj-x) = g'j-i • (cj + oy-i ~~ x ) an d 
/(ay-i) - f{x - Cj) < g'j_ x \cj + ay_i - x); then f * gj{x) - f * gj-i(x) > 0. 

The second statement can be proved similarly. □ 

Another formulation of Lemma 4.4 is that gj > f*gj-i and that > f*gj 
and that the two fonctions intersect at least once. Hence gj and cannot 
appear in the minimum of f *gj and / * gj—i- By transitivity, there is no need 
to compute entirely the convolution of the convex function by every affine 
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Figure 2. Convolution of a convex function by a concave function. 



component of the decomposition of the concave function. If there are more 
than two segments, successive applications of this lemma show that only the 
position of the segments of the concave function must be computed, except 
for the extremal segments. 

The following lemma shows that / * gj and / * 9j-i intersect in one and 
only one connected component, as for a given abscissa, the slope of f * gj is 
less than the one of / * gj-i- 

Lemma 4-5. — Let i£R where f * gj and f * gj-i are defined and differ- 
entiate. Then 

( 41 ) -j^/ *Si(*)< ^f* gj -i(x). 

Proof. — For x G [o , ay], / * gj(x + cy) = f(x - Cj) + g(cj) and / * gj-x(x + 
Cj-i) = f(x) + g(cj-i). As / is convex and Cj > Cj—x, the result holds on 
[ao + Cj,otj + Cj]. Similarily, the result holds for x G [a.j-\ + Cj, a n + Cj]. 

On [cj + atj, Cj + otj-i], f * gj is composed of segment gj concatenated with 
the segments /j, i G [atj-\,oij], possibly troncated on the right and / * gj-\ 
is composed of segments fi, i G [aj—i,atj] concatenated with gj-i, possibly 
troncated on the left. As Vi G [ay_i, ctj], g'j_ 1 < f- < g'p the result holds. □ 

If one sets by convention ^/ * gj(x) = —oo for x < Cj_i + ao and 4^/ * 
gj{x) = +oo for x > Cj + a n , then the inequality always holds. 

The intersection of / * gj-i and / * gj can then happen in one and only one 
of the four cases: 

1. g)_ x and g c - intersect; 3. and g^ intersect; 

2. gj_ x and g? intersect; 4. g^_ x and o| intersect. 
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The following theorem is another consequence of these lemmas and is more 
precise about the shape on the convolution of a convex function by a concave 
function. 

Theorem 4-6 (convolution of a convex function by a concave func- 
tion) 

The (min, plus) ^-convolution of a convex function by a concave function can 
be decomposed in three (possibly trivial) parts: a convex function, a concave 
function and a convex function. 

Proof. — We use the notations defined in the former lemmas. 

We now show by induction that the convolution of / by minfc<j gj , denoted 
hj, is composed of 

(i) a convex part hj, which is the a restriction of go(x) = f(x — cq) + g(co) 
to [co + a , f3j], with (3j < c + a n ; 

(ii) a concave part hj, which is a minimum of some segments g^, k < j (up 
to some translation) defined on [/3j,7j]; 

(iii) a convex part h?, gj{x) = f(x - c j+ i) + g(c j+ i) for x G [lj,a n + Cj+i], 
7j > «o + Cj+i. 

Note that with these conventions, the reals (3j and 7j are uniquely deter- 
mined at each step of the induction. 

The case with j = is a direct consequence of Lemma 4.1. The case with 
j = 1 is a consequence of Lemmas 4.4 and 4.1. The graphs of / * g\ and f * go 
intersect once and only once (where they are defined), and in [c\ + a.\,ci + ao\- 
Depending on when this intersection occurs, the concave part will be trivial, be 
made of only one (part of a) segment of g, or a minimum of the two segments 
go and g\. 

Suppose now that the result holds for hj and consider hj+i = min(/ij,/ * 
gj+i). The argument is exactly the same as for j = 1: hj and f*gj+% can only 
intersect once and only once. Indeed, hj is the minimum of functions such 
that 35/ * Sk{x) > ^f*g j+ i(x), and then j^hj(x) > £,f * g j+l (x). Note 
that this intersection has to occur after the point Cj+i + ay+i. 

Moreover, as gj does not intersect <?J +1 and that hj is a part of g 2 j (by the 
induction hypothesis), hj does not intersect 9j + i- 

Therefore, only one of the four following cases may occur. 

1. h c j intersects g c j+l and h l j+1 = hj, h c j+1 = mm(h°j,g] +1 ), h 2 j+1 = g 2 j+l , 
f3j + i = Pj and 7j+i = a j+1 + c j+2 - 

2. h°j intersects gj +1 at y and h l j+1 = h), h c j+1 = h C j, h 2 j+l = g 2 +1 , = f3j 
and 7j + i = y. 

3. h) intersects g c j+l at y and hj +1 = hj, h c j+1 = h C j, h 2 +l = g 2 +1 , (3 j+1 = y 
and 7j+i = a j+ i + c j+2 . 
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4. hj intersects at y and = hj, h c - +l is trivial and = 

□ 



If the concave function is composed of m segments and the convex function 
of n segments, then the convolution of those two functions can be computed 
in time 0(n + mlogm). The logm term comes from the fact that one has 
to compute the minimum of m segments (see [BT08] for more details). If 
the functions are now defined on N, then, as no intersection point has to be 
computed for the minimum, the time complexity is 0{n + m). The correspond- 
ing algorithm is given in Algorithm 2, where without loss of generality (the 
(min,plus)-convolution is shift-invariant), the functions / and g are defined on 
N and finite between respectively and n, and and m. The slopes of the 
functions are thus f[ = f(i) — f(i — 1) and g\ = g(i) — g(i — 1). 



Algorithm 2: Convolution of a convex function by a concave function 
Data: / : [0, n] — > ffi a convex function with slopes (r^), g : [0, m] — > M. a 

concave function with slopes (pi). 
Result: h = f * g 
begin 

/* Initialization */ 

k <- 0; 

while ^m + ndo h(k) < — hoo; k <— k + 1; 
i^0;j^0; fc(0)«-/(0)+s(0); 

/* First convex part of the convolution */ 
while fl ^ g' do 

ii^i + l; h(i) <-/(*) + 5(0); 

/* Concave part of the convolution */ 

j<-j + l\ h(i+j) <- f(i)+g(j); 
while j < m do 

while g'- < f' i _ 1 do i ^— i — 1; 

h(i + j) <r- mm(h(i + j),f(i)+g(j)); 

h{i + j + 1) i- mm(h(i + j + 1), /(») + g(J + 1)); 

L i <- 3 + 1; 

/* Second convex part of the convolution */ 
while i < n do 

i -<— i + 1; /i(i + m) ^— min(/i(i + m), /(i) + g(m)); 



FAST WEAK-KAM INTEGRATORS 



25 



4.2. Application. — We now go back to our initial problem. As already 
explained above, the computation of T t T £ u(t, x) given in (32) is made of two 
steps: 

— (min,plus)-convolution of u and h : x i— )■ tK*(*); 

— subtract rV(t,x). 

Note that the (min,plus)-convolution described above is here defined on 
functions that have a non-bounded support. But in the periodic case, u is 
1-periodic and h is convex with a global minimum. Then, to compute the 
convolution, it is enough to compute it on a single period (the (min,plus) 
convolution preserves the periodicity), and replace h by its restriction on a 
support of size 2 centered on its minimum. If e = 1/N with the notation of 
the previous section, then both functions u and h are defined on grids of size 
iV and 2N respectively. 

The convolution of u and h can be efficiently computed following these steps: 

1. Decompose u into convex and concave parts. This can be done in linear 
time: the three first points determine if a part is concave or convex. 
Then, this part is extended as much as possible while preserving the 
concavity or convexity and so on. 

2. For each convex or concave part, perform the convolution with h using 
Algorithms 1 or 2. 

3. Take the minimum of all these convolutions. 

The complexity of this Algorithm is then 0(cN), where c is the number of 
components in the decomposition of u into concave/convex parts. 



4.3. Implementation issues. — The main issue with this algorithm is that 
c - the number of components in the decomposition of u - can become very 
large, and then lead to a quadratic time complexity, which is the complexity 
of a naive algorithm for computing the convolution. Experimentally, the rea- 
son for this is that, due to the discretization of u, nearly affine parts, after 
performing the convolution several times, are computed as fast alternations of 
convex and concave parts. As shown in Figure 3, one solution to make the 
computations more efficient would be to consider those parts are convex and 
use Algorithm 1. 

To do this, we decompose u into convex and concave parts with a tolerance 
(we do not request for convex parts to have increasing increments, but the 
increments to have an increase more than —rj). We will discuss this in the 
next section. The choice of an optimal tolerance rj, as well as the comparison 
with parallel implementations, will be the subject of further studies. 
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Figure 3. Approximation of the convolution: plain line shows the 
convolution of u and h, and the dashed line shows the function com- 
puted using Algorithm 1 when u is not convex, but has very small 
variations. 



5. Numerical simulations 

We take the Hamiltonian (3) with P = 1 and V(t, x) = 1 — cos(x) on [—it, it]. 
We take N = 600 grid points, that is e = 1.7e — 3, and r = \fe = 0.04. 
The comparisons are made with the fifth-order WENO algorithm (WEN05 
see [JS96]) with 100 grid points, which will be considered here as the exact 
solution. 

In a first simulation, we calculate the solution at times t = 1,2,5 and 15 
with initial value u(0,x) = cos(2x). We see the very good agreement between 
the solution given by the WEN05 algorithm. However, the CPU time required 
by our algorithm is about 3 times the CPU time of the WEN05 algorithm 
(of order 2s to reach t = 15). In this case, for 600 grid points in [—it, it], 
the number of convex/concave components c in the decomposition of u is of 
order 100 at each time step. The plot, rescaled such that u(t, —it) = 0, is 
shown in Figure 4. Note that after the time t = 20, the solution has converged 
and remains constant for larger times. According to Remark 3.8, the solution 
observed is thus very close to the weak-KAM solution u*(x). 

In a second simulation, we use the approximated convolution with a tol- 
erance r\ = 10 x e 2 = 2.7e — 5. In this case, the number of convex/concave 
components c is always of order 10 (and equal to 3 - as expected from the 
shape of the solution - when the stationary state is attained), and the CPU 
time is reduced by a factor 20 when compared with r] = (and about 4 times 
quicker than the WEN05 algorithm), without any significant deterioration of 
the accuracy. We show the result in Figure 5. 

Finally, we take P = 2 and the potential function V(t,x) = sin(t) cos(2x), 
e = 0.01 and r = 0.1. We take r\ = 0, and consider the initial data uq{x) = 
— cos(3x). In Figure 6 we plot the evolution of u(t, —it) with respect to the 
time. We observe the linear growth predicted by the result of the previous 
section, for a time interval [0, 10000]. Note that there are oscillations in the 
evolution of u(t, —it), but at a scale too small to be visible on the plot. 
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t=1 t=2 




Figure 4. Solutions with n = 0. CPU = 6s at t = 15 

Appendix: Proof of the a priori compactness Proposition 3.2 

We start with a lemma. 

Lemma 5.1. — Assume that the hypothesis (i) and (Hi) are satisfied. Recall 
that L(t,x,v) = K*(v) — V(t,x). For any T > 1, there exists a constant V 
such that for any x, y G W 1 and T > and t > 1, if \x — y\/t < T and 7 that 
minimizes the quantity 

inf / L(T + s,7(s),7(s))ds, 

7(0)=x Jo 

i(t)=y 

then 

VO^a^a + l^t, | 7 (a) -7(a+ 1)| < T'. 

Proof. — Without loss of generality, we will assume that L is positive. Indeed, 
the potential V is bounded, and adding a constant doesn't change the mini- 
mizers. In this case, and under the hypothesis (i), there exists a nonnegative, 
increasing function a which tends to 00 at 00, such that 

Vt,x,v, L(t,x,v) ^ a(\v\)\v\. 
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Figure 6. Evolution of u(t, —ir) over long times 



The idea of the proof is that if 7 at some point has a great velocity, then it 
must be slow later. It is then better to "slow down" the fast part and accelerate 
the "slow" one. 

First, we set some notations. For all (x,y,t,T) £ 1" x R" x R + x R, let 
A t T (x,y)= inf / L(r + s, 7 ( S ),7( S ))d S , 

7(0)=x Jo 

i(t)=y 

be the Lagrangian action. 
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We start by showing, that the superlinearity of L implies a superlinearity 
of Asp. As already done a few times, we may bound the action by comparing 
with a straight line, using that L is uniformly bounded on sets of the form 
Kxl"x B(0, R), where B(0, R) is the ball of radius R > in R n : 

, , f* , (t — s)x + sy y — x^, 

Atp(x,y) < / L{T + s, K - >- )ds 

Jo it 

(42) < max (i^i, l), 

for some increasing function z i— > C + (z) defined on IR + . Let 7 realizing the 
infimum, and set 

£ = {sG[0,i] such that |'y(a)| > ^^ }- 
Then we get, using that L > 0: 

(43) ^(ar,y)= t L(T + s,j(s),j(s))ds 

Jo 

> ^L(T + S ,7( S ),7( S ))ds 

'i — yk /" . /I^ - y\\\ x ~ y\ 



2i ' J e v 2t 7 2 

The last inequality comes from the fact that when 7 is going at speed less 
than \x — y\/2t for time t, it cannot travel more than \x — y\/2, therefore the 
integral is greater than \x — y\/2 by the triangular inequality. In other terms, 
equations (42) and (43) state that there are two positive functions C + and 
which can be easily made increasing, such that 



(44) C -i}l^A)\l^A ^ A T^V) ^ C+(i^)max(^,l 



it 

hp 

t ' t t " ' v t ' V t 

Moreover, thanks to the superlinearity of L those functions go to 00 at 00. 
Now consider 

- T" > r such that 20C+(r) < c-(T"), 

- V" > T" such that T"'/T" € N* and 30C+(20r") < C-(T"'), 

- and finally V > V" such that 40T'" /T" < C~ (T')/C + (T). 

Let us verify that V satisfies the requirements of our lemma. 

Assume by contradiction that for some x,y £ M n , t, T G M + such that 
\x — y\ < tT and 7 realizing the action Atp(x, y), there is an a £ [0, t — 1] such 
that (7(a) — j(a + 1)| r'. As 7 is a minimizer, we have (using that L > 



30 



ANNE BOUILLARD, ERWAN FAOU & MAXIME ZAVIDOVIQUE 



and r > 1) 

r-a+l 



(45) tTC+(T)^A^{x,y)^ f L(T + s, 7 (s),j(s))ds 

J a 

^Ai +a ( 7 (a), 7 (a + i)) ^r'c-(r'). 

Hence we obtain 

rc-rr) 4or w 

( 46 ) * > ^ttW > 



rc+(r) r" ' 

using the fact that r' > T and the definition of V . 

We now assume that a < t/2, the other case may be treated similarly. Let 
b £ [a, a + 1] be the smallest such that \"f(a) — 7(6) | = T'", and consider the 
sequence 

Ci = 6 + 2f — , i 6 {(),••• ,A}, 

where A; is greatest possible integer such that ^ t. Note that using (46) 
and a ^ t/2, we have A; 9, and that for i S {0, • • • , A}, we have q+i — Cj = 
2T'"/T". 

We claim that there exists an iq £ {0, • • • , A;} such that 

l7(ci ) -7(ci +i)| 



< r". 



Indeed, otherwise we would have, using (44) 

^4(x,y) r +1 L(T + S ,7( S ), 7 ( S ))d S 

i=0 Ci 

^ \p l7(c»o) ~ 7(c»o+i)l c -( r //) 
i=0 c *o+l ~~ c «o 

> E^r^n- 

i=0 

By definition of k, we have that (A; + 1) x 2T"'/T" > t/2 - 1 while using 
(46), we have t > 40r /// /T // ^ 40 and 2r /// /r // < t/18. Hence we deduce that 
k x 2T"'/T" > i/3. As r" ^ T, the previous equation yields 

Air(x,y) > trC-(T") > tTC + (F), 

which is absurd in view of (45). 

Now we find a contradiction by constructing a curve 5 which has an action 
less than 7. Let [c,d] = [cj ,Cj +i]. Recall that N := Y"'/T" is an integer. We 
define the curve 5 as follows: 
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— 5(s) = 7(s) if s € [0, a] U [d, t]; 

— on [a, b + N], 5 coincides with the curve minimizing A b ^!~^ l a ~ a (7(a), 7(b)) ; 

— on [b + N, c + N], 5 is the translate of 7 : 5(s) = 7(5 — AT); 

— on [c+N, d] (recall that d = c+2N) 5 coincides with the curve minimizing 
A» +c+N (^c)Md)). 

We now compute the difference of action between 7 and 5, recalling that L is 
1-periodic in time: 

I L(T + s,j(s),j(s))ds- f L(T + s,5(s),5(s))ds 
Jo Jo 

rb rd 

= / L(r + s,7( S ), 7 ( S ))d S + / L(T + s, 7 (s),j(s))ds 

J a J c 

/•b+N pd 

- / L(T + s,5(s),6(s))ds- / L(T + s, S(s), 8{s))ds 

Ja Jc+N 

^T"'c-(v'") + ^-T"c-(r") 

-pill -pill 

- ^r"c+(r") - — (2r")c + (2T") > o. 

This contradicts the minimality of 7. 

□ 



Remark 5.2. — In the previous proof, we only used the fact that L is peri- 
odic in time. In [Itu96] , a similar result is proved when L is periodic in space 
(instead of in time). The idea of the proof is the same except that, when 
constructing the curve 8, instead of translating it in time (in third part of the 
construction), it is translated in space, while the "fast" part of 7 between a 
and b is replaced by a geodesic (straight line) between 7(a) and the closest 
point from 7(a) in the grid 7(6) + Z n . 

We now prove the lemma 3.2: 

proof of lemma 3. 2. — Recall now that L is periodic both in time and in space 
and that its Euler-Lagrange flow is complete. As in the previous lemma, 
assume L > 0. Let V and V be as in the previous lemma, and 7 be a minimizer 
such that (7(0) — 7(i)|/i ^ T. The curve 7 is then a trajectory of the Euler- 
Lagrange flow. Let moreover O^a^a + l^t. Finally, by superlinearity of 
L, let A{1) be given by Equation (8), such that L(t,x,v) ^ \v\ — A(l). We 
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therefore obtain that, with the notations used in the previous proof, 



By periodicity of the Lagrangian, and completeness of the Euler-Lagrange 
flow, there exists a constant D' depending only on D, such that I7I ^ D' on 
[a + sq — 1, a + So + 1] l~l [0, t] D [a, a + 1]. Since a is arbitrary, this finishes the 
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